home *** CD-ROM | disk | FTP | other *** search
/ Games of Daze / Infomagic - Games of Daze (Summer 1995) (Disc 1 of 2).iso / x2ftp / msdos / math / newmat08 / newmat1.cpp < prev    next >
C/C++ Source or Header  |  1995-01-11  |  4KB  |  140 lines
  1. //$$ newmat1.cpp   Matrix type functions
  2.  
  3. // Copyright (C) 1991,2,3,4: R B Davies
  4.  
  5.  
  6. #include "include.h"
  7.  
  8. #include "newmat.h"
  9.  
  10.  
  11.  
  12. /************************* MatrixType functions *****************************/
  13.  
  14.  
  15. // all operations to return MatrixTypes which correspond to valid types
  16. // of matrices.
  17. // Eg: if it has the Diagonal attribute, then it must also have
  18. // Upper, Lower, Band and Symmetric.
  19.  
  20.  
  21. MatrixType MatrixType::operator*(const MatrixType& mt) const
  22. {
  23.    int a = attribute & mt.attribute & ~Symmetric;
  24.    a |= (a & Diagonal) * 31;                   // recognise diagonal
  25.    return MatrixType(a);
  26. }
  27.  
  28. MatrixType MatrixType::SP(const MatrixType& mt) const
  29. // elementwise product
  30. // Lower, Upper, Diag, Band if only one is
  31. // Symmetric, Valid (and Real) if both are
  32. // Need to include Lower & Upper => Diagonal
  33. // Will need to include both Skew => Symmetric
  34. {
  35.    int a = ((attribute | mt.attribute) & ~(Symmetric + Valid))
  36.       | (attribute & mt.attribute);
  37.    if ((a & Lower) != 0  &&  (a & Upper) != 0) a |= Diagonal;
  38.    a |= (a & Diagonal) * 31;                   // recognise diagonal
  39.    return MatrixType(a);
  40. }
  41.  
  42. MatrixType MatrixType::i() const               // type of inverse
  43. {
  44.    int a = attribute & ~(Band+LUDeco);
  45.    a |= (a & Diagonal) * 31;                   // recognise diagonal
  46.    return MatrixType(a);
  47. }
  48.  
  49. MatrixType MatrixType::t() const
  50. // swap lower and upper attributes
  51. // assume Upper is in bit above Lower
  52. {
  53.    int a = attribute;
  54.    a ^= (((a >> 1) ^ a) & Lower) * 3;
  55.    return MatrixType(a);
  56. }
  57.  
  58. MatrixType MatrixType::MultRHS() const
  59. {
  60. // reomve symmetric attribute unless diagonal
  61.    return (attribute == Dg) ? Dg : (attribute & ~Symmetric);
  62. }
  63.  
  64. Boolean Rectangular(MatrixType a, MatrixType b, MatrixType c)
  65. {
  66.    return
  67.       ((a.attribute | b.attribute | c.attribute) & ~MatrixType::Valid) == 0;
  68. }
  69.  
  70. char* MatrixType::Value() const
  71. {
  72. // make a string with the name of matrix with the given attributes
  73.    switch (attribute)
  74.    {
  75.    case Valid:                              return "Rect ";
  76.    case Valid+Symmetric:                    return "Sym  ";
  77.    case Valid+Band:                         return "Band ";
  78.    case Valid+Symmetric+Band:               return "SmBnd";
  79.    case Valid+Upper:                        return "UT   ";
  80.    case Valid+Diagonal+Symmetric+Band+Upper+Lower:
  81.                                             return "Diag ";
  82.    case Valid+Band+Upper:                   return "UpBnd";
  83.    case Valid+Lower:                        return "LT   ";
  84.    case Valid+Band+Lower:                   return "LwBnd";
  85.    default:
  86.       if (!(attribute & Valid))             return "UnSp ";
  87.       if (attribute & LUDeco)
  88.          return (attribute & Band) ?     "BndLU" : "Crout";
  89.                                             return "?????";
  90.    }
  91. }
  92.  
  93.  
  94. GeneralMatrix* MatrixType::New(int nr, int nc, BaseMatrix* bm) const
  95. {
  96. // make a new matrix with the given attributes
  97.  
  98.    Tracer tr("New"); GeneralMatrix* gm;
  99.    switch (attribute)
  100.    {
  101.    case Valid:
  102.       if (nc==1) { gm = new ColumnVector(nr); break; }
  103.       if (nr==1) { gm = new RowVector(nc); break; }
  104.       gm = new Matrix(nr, nc); break;
  105.  
  106.    case Valid+Symmetric:
  107.       gm = new SymmetricMatrix(nr); break;
  108.  
  109.    case Valid+Band:
  110.       {
  111.          MatrixBandWidth bw = bm->BandWidth();
  112.          gm = new BandMatrix(nr,bw.lower,bw.upper); break;
  113.       }
  114.  
  115.    case Valid+Symmetric+Band:
  116.       gm = new SymmetricBandMatrix(nr,bm->BandWidth().lower); break;
  117.  
  118.    case Valid+Upper:
  119.       gm = new UpperTriangularMatrix(nr); break;
  120.  
  121.    case Valid+Diagonal+Symmetric+Band+Upper+Lower:
  122.       gm = new DiagonalMatrix(nr); break;
  123.  
  124.    case Valid+Band+Upper:
  125.       gm = new UpperBandMatrix(nr,bm->BandWidth().upper); break;
  126.  
  127.    case Valid+Lower:
  128.       gm = new LowerTriangularMatrix(nr); break;
  129.  
  130.    case Valid+Band+Lower:
  131.       gm = new LowerBandMatrix(nr,bm->BandWidth().lower); break;
  132.  
  133.    default:
  134.       Throw(ProgramException("Invalid matrix type"));
  135.    }
  136.    
  137.    MatrixErrorNoSpace(gm); gm->Protect(); return gm;
  138. }
  139.  
  140.